Scaling behavior in a stochastic self-gravitating system 
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Abstract 

A system of stochastic differential equations for the velocity and density of a classical self- 
gravitating matter is investigated by means of the field theoretic renormalization group. The exis- 
tence of two types of large-scale scaling behavior, associated to physically admissible fixed points of 
the renormalization-group equations, is established. Their regions of stability are identified and the 
corresponding scaling dimensions are calculated in the one-loop approximation (first order of the 
e expansion). The velocity and density fields have independent scaling dimensions. Our analysis 
supports the importance of the rotational (nonpotential) components of the velocity field in the 
formation of those scaling laws. PACS numbers: 05.45.-a; 05.10.Cc; 04.40.-b. 
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The present Universe at short and moderate distances 
is inhomogeneous, being filled by numerous structures 
of many scales, from galaxies to galaxy clusters and su- 
perclusters. On the contrary, at larger scales or earlier 
stages the Universe is generally taken nearly homoge- 
neous and isotropic. It is believed that gravity reinforces 
small asymmetries in the velocity and density fields, and 
the structures observed today are due to instabilities in 
an initially uniform self-gravitating medium. The first 
structures to form are "pancakes," thin in one dimension 
and of large extent in the two others. Further evolution 
and interaction of the pancakes develops patterns with 
complex (fractal or honeycomb) geometry in the distribu- 
tion of matter. Excellent reviews of classical cosmology 
are given in Refs. [1-3]. 

The full relativistic treatment of the large-scale struc- 
ture formation, especially of its nonlinear stage, is an 
extremely difficult task. Therefore, the development of 
instabilities is usually studied within the framework of 
simplified dynamical models of a classical self-gravitating 
fluid (or system of particles): Vlasov-Poisson model, ad- 
hesion model and its modifications with different types 
of nonlinearity, pressure and viscous terms and random 
forces, Boltzmann equation (or A-body simulations of 
dark matter) and so on; see e.g. Refs. [3,4] for reviews 
and discussion. 

The link between the complex geometrical structures 
and the coarse-grained (hydrodynamic) description is 



provided by nontrivial scaling behavior exhibited by cor- 
relation functions of the density or velocity fields, such as 
the galaxy-galaxy correlation function £(r); see e.g. [5]. 
At large scales it reveals a power-law behavior £(r) oc 
r ' , where 7 is determined from catalogs to be between 
1.3 and 2.1 for r of order of the Megaparsec; see e.g. [6-10] 
and references therein. The most recent studies give the 
values between 1.6 and 1.9, in particular, 7 = 1.75 ±0.03 
according to [9]. 

The scope of theory is to derive such behavior on 
the basis of an appropriate dynamical model, to inves- 
tigate the universality of the exponent 7, that is, its 
(in) dependence on the model parameters, and to calcu- 
late 7 within a consistent approximation or systematic 
perturbation scheme. 

Scaling laws are typical of equilibrium phase tran- 
sitions, and the most adequate tool to study them is 
the renormalization group (RG). It is also applicable 
to nonequilibrium dynamical phenomena as disparate as 
surface growth, random walks, nonlinear diffusion and 
turbulence; see e.g. [11] for a review. For a given model, 
the RG allows one to prove the existence of scaling 
regime(s), to determine the range of its stability in the 
space of model parameters, and to calculate the scaling 
dimensions in the form of regular perturbation series (e 
expansions). 

The RG approach to the problem of self-gravitating 
medium was pioneered in [12-14]. In these studies, the 
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full set of equations (hydrodynamic equation for the New- 
tonian fluid, continuity equation and the Poisson equa- 
tion for the gravitation force, the system known as the 
Vlasov-Poisson equations) was reduced to a single equa- 
tion for a purely potential velocity field. The resulting 
equation (similar to the well-known stochastic Burgers 
equation but with a time-dependent "mass" term) was 
augmented by a Gaussian random force (noise) that rep- 
resents the influence of fluctuations and dissipative pro- 
cesses on the evolution of fluid, arising from viscosity, tur- 
bulence, explosions, gravitational waves and so on. The 
dynamical RG approach of [15-17] was then adopted to 
derive the scaling regimes and exponents. With an appro- 
priate choice of the parameters, the model reveals scal- 
ing behavior with a nonuniversal exponent 7; its value 
depends on the characteristics of the forcing and can be 
adjusted to the value 7 « 1.7 [14] in agreement with the 
observations. 

In spite of this obvious success, the analysis of Refs. 
[13,14] raises serious questions about its internal con- 
sistency and interpretation of the results. It is well 
known that the stochastic Burgers (or Kardar-Parisi- 
Zhang) equation has no infrared-attractive fixed point 
in the physical range of parameters within the e expan- 
sion. This fact immediately follows from the first-order 
expressions of Refs. [15-17]. It was confirmed by the 
two-loop calculation [18] and then proved to all orders of 
the perturbation theory [19]. The existence of a strong- 
coupling fixed point, although supported by numerical 
simulations, remains an unproved hypothesis. 

The authors of [14] studied an extended ("massive") 
version of the KPZ model, and the only attractive fixed 
point revealed in their analysis corresponds to nonzero 
value of the "mass." In fact, its value is comparable with 
the largest, ultraviolet, momentum scale of the problem, 



so that the mass term at the fixed point is by no means 
small. Therefore, in the spirit of the Landau theory, the 
viscous and nonlinear terms in the equation should be 
discarded as infrared-irrelevant. This leads to a simple 
Gaussian model whose critical exponents are easily found 
exactly. They agree with the answers obtained in [14] 
for the full (interacting) model, but the situation on the 
whole is not quite satisfactory. It looks rather strange 
that the nonlinearity does not play any role in shaping 
the large-scale asymptotic form of the correlation func- 
tions, and that the latter are so directly determined by 
the choice of the random forcing. The absence of the 
appropriate fixed point for the massless (critical) model 
suggests that such results can be very sensitive to the ap- 
proximations made in the derivation of the model. The 
most important of them is the assumption of parallelism, 
which reduces the system of equations for the density 
and velocity fields to an equation for a single, purely po- 
tential, velocity field. More extensive analysis of the full 
problem would therefore be desirable. 

The original (deterministic) set of Vlasov-Poisson 
equations for the velocity Ui(t, x), density p(t, x) and the 
gravitational potential tp(t, x) in the comoving frame of 
reference reads (see e.g. [3,4]): 

adtUi + aui + (ujdj)ui — diip = 0, 
ad t p + di(pui) = 0, d 2 ijj = 4nGa 2 (p - p ), (1) 

where G is the gravitational constant, the cosmic scale 
factor a(t) is a prescribed function of the proper time, p n 
is the mean density and d 2 = didi is the Laplace opera- 
tor. 

We change to the new variables fa = Ui/a, 9 = 
c 2 (p — po)/po an d c 2 = AitGa 2 po and rewrite the system 
(1) in the form: 



(d t + fadj) fa = -Hfa + v (S i:j d 2 - didj) fa + u Q v didj4>j - cdid 2 9 + f it d t 6 = v is d 2 9 - di(fa9) - c(d l fa). (2) 



Here, H = a/ a is the Hubble function and d~ 2 is the 
Green function of the Laplace operator. We have elimi- 
nated the potential tp using the last equation in (1), as- 
sumed that a/a -C d t u/u, and added viscous terms and 
the random force fi(t, x). We stress that the velocity field 
fa oc Ui is not purely potential, so that two independent 
viscosity coefficients vq and uqi/q have been introduced 
in the equation for fa. The viscous term VQV d 2 9 is usu- 
ally not included in the continuity equation, but it is not 
forbidden by dimensional reasons or symmetry and thus 
is needed to ensure the renormalizability of the model. 
Then the RG equations should be solved with the physi- 
cal initial condition v — 0, but if the IR attractive fixed 
point is unique, the large-scale behavior will be the same 
as for nonzero Vo; cf. the discussion in [20]. Dimen- 
sional analysis shows that the pressure term is infrared- 
irrelevant (in the sense of Wilson) in comparison to the 



gravitational force and thus it was dropped in (2). 

The hydrodynamic description of the properly 
smoothed (coarse-grained) fields and, in particular, in- 
clusion of the viscous terms and random forcing, can be 
justified by various arguments [21-23]. For simplicity, we 
shall neglect the time dependence of the viscosity coeffi- 
cients, suggested by those studies, treating it as a kind of 
second-order effect: the viscosity coefficients are "small" 
and their time dependence is "slow." In contrast to equi- 
librium systems, there is no universal relation between 
the viscosity coefficients and the correlation functions of 
the random force. We shall take it Gaussian, white in 
time (this is necessary to ensure the Galilean symme- 
try of the stochastic problem (2)), with zero mean and a 
given correlator 

(Mt,\L)fitf, -k)) = S(t - t')D(k) {P tJ + aQ y } . (3) 
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Here Pij = Sij —kikj/k 2 and Qij = kikj /k 2 are the trans- 
verse and the longitudinal projectors, a > is an arbi- 
trary parameter and D(k) is a function of the modulus of 
the wave vector k = |k|. The simplest possible choice is 
D(k) = Z? =const (spatially decorrelated forcing). An- 
other possibility, widely used in models of nonequilibrium 
critical phenomena, is a power-law correlation function: 
D(k) = D' k i - d - 2 ^- see e.g. [15,17,24-27]. Here d is the 
(arbitrary) dimensionality of space and r\ an arbitrary 
exponent; the notation is explained by convenience rea- 
sons. In what follows, we shall refer to these two cases as 
the local and nonlocal ones. They are closely related for 
d < 4 (see below), but it is instructive to discuss them 
separately in the beginning. 

The RG analysis of a stochastic problem like (2), (3) 



where Df is the correlator (3) and the needed integra- 
tions over t, x and summations over the vector indices 
are implied. 

The field theoretic formulation means that the corre- 
lation functions of the stochastic problem (2), (3) can 
be represented as functional averages with the weight 
exp5'($) with action (4). This allows one to use a well- 
developed formalism (power counting plus symmetries of 
the model) to analyze the relation between the IR and 
UV problems and the UV renormalizability of the model. 
For the local case, it shows that the upper critical dimen- 
sion for the model is d — 4: the nonlinearity in (2) is IR 
irrelevant for d > 4 (perturbation theory is applicable, 
no scaling and universality are expected). For d < 4, the 
terms of the ordinary perturbation theory suffer from IR 
singularities and cannot be used to describe the large- 
scale behavior of the problem. For small e = d — 4, 
the problem of the IR singularities is closely related to 
that of the UV divergences (poles in e). The latter is 
solved by the standard UV renormalization procedure: 
it shows that the model (4) is multiplicatively renormal- 
izable, that is, all the poles in s in its correlation func- 
tions are removed by the rescaling of the fields <f> and the 
parameters Dq, vq, uq, vo, a (the proof of this statement 
is the most nontrivial stage of the analysis). The arbi- 
trariness in the renormalization procedure leads to the 
RG equations: first-order differential equations for the 
correlation functions with coefficients calculated within 
the ordinary perturbation theory. In order to draw any 
definite conclusions from the RG equations, one has to 
calculate their coefficients at least in the simplest (one- 
loop) approximation. We performed the calculation and 
found out that, in contrast to the Burgers or KPZ mod- 
els, the RG equations of the extended model (4) have 
the only IR attractive fixed point in the physical range 



includes four important steps: field theoretic formula- 
tion; analysis of its renormalizability; derivation of the 
corresponding RG equations; analysis of the fixed points 
of these equations. This analysis for our problem is tech- 
nically involved (already in the simplest one-loop approx- 
imation) and will be presented elsewhere, along with the 
details of the practical calculation. In many respects, it 
is close to the field theoretic RG analysis of the stochastic 
Navier-Stokes equation [24-27] and especially to the case 
of a strongly compressible fluid studied in [20] . Below we 
only give the main points and conclusions. 

According to a general theorem (sec e.g. [11,26]), the 
stochastic problem (2) is equivalent to the field theoretic 
model of a doubled set of fields <f> = {</>', 6', <fi, 0} with 
action functional: 



(4) 

of parameters (the ratios of the viscosity coefficients and 
the amplitude factors in pair correlation functions are 
positive) . 

This means, in particular, that in the IR range (the 
scales large in comparison to the typical UV length scale, 
built of D and v 0l and times large in comparison to the 
corresponding time scale), the correlation functions of 
the velocity <fi and the density (more precisely, of the 
field 9 = c 2 (p — po)/po have a scaling (self-similar) form: 

(0(t, x)0(t + r, x + r)> ~ r- 2A * F4. . .) , 
(0(t,x)0(t + T,x + r)) ~r- 2A *F e (...), (5) 

where r = |r| and the scaling functions F v g depend on 
(critically) dimensionless variables t-t At , H-r A " , c-r Ac 
(for the equal-time correlation functions, the first vari- 
able is absent). The dimensions A are universal in the 
sense that they are independent of the values of the pa- 
rameters uo,vo,a and can be calculated as scries in e. 
The first-order (one-loop) calculation gives: 

A = 1 - e/2, A e = A c = 2 - e/2, 

A T = -2 + e/2, A ff = 2 + e/2, (6) 

with corrections of order e 2 and higher. From represen- 
tations (5) it follows that under the rescaling 

r^r/A, t^tA a -, H -> HA Ah , c -» cA A " (7) 

with arbitrary A > 0, the correlation functions behave as 

(w) -» A 2A " (w) , (00) -» A 2Ae (90) . (8) 

The formulation (7), (8) is in fact more general because 
it remains true if the parameters H, c depend on t, x (in 
the original problem they indeed depend on t), while the 
more explicit formulae (5) imply that they are treated as 
constants. 



S($) = {\/2)4fDj4f + <(>'.{-(% + fadj +H)<f>i + is (M 2 - didj) <f>j + uovodidjfa - cd l d- 2 9) 
+ 9' {-8 t 9 + v u d 2 9 - di(<f>i0) - c{di<Pi)} , 
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The RG representations (5) are the result of certain in- 
finite resummation of the primitive perturbation theory, 
that is, of the expansion in the nonlinearity in Eqs. (2) 
around the zero-order (Gaussian) approximation. In our 
case, however, the latter is unstable with respect to any 
small perturbation, as is easily seen from the fact that, for 
c 2 > 0, the retarded zero-order response function grows 
in time and, as a result, perturbative diagrams contain 
infrared divergences. However, it can be argued that this 
instability does not hinder the use of the RG in studying 
the self-similar behavior. The parameters H and c in (4) 
have integer positive dimensions and, in this respect, they 
are analogous to masses (in the language of the quantum 
field theory) or to the deviation of the temperature of 
its critical value, AT = T — T c (in models of critical be- 
havior). From the general theory of UV renormalization, 
it is well known that the UV divergent parts of the dia- 
grams are polynomials in such "IR relevant parameters." 
Therefore, they can be calculated for AT > (or c 2 < 
in our case), where the terms of the perturbation theory 
are finite, and then extrapolated to the region AT < 
(or c 2 > 0). There, the ordinary perturbation expan- 
sion ceases to make sense due to IR divergences and one 
should either change to an improved perturbation the- 
ory (for critical behavior) or to consider a nonstationary 
problem (for a self-gravitating system). It is important 
here that this rearrangement does not affect the UV di- 
vergent parts of the correlation functions (and hence the 
counterterms) . These arguments show that the critical 
exponents (calculated from the UV counterterms) are the 
same below and above T c , and, in our case, they support 
the scaling relations (5) and (7), (8) with the dimensions 
(6) for the "unstable" case c 2 > 0. 

For the nonlocal case, that is, D(k) oc k 4 ~ d ~ 2r > in (3), 
analysis shows that the model (4) appears multiplica- 
tively renormalizable if d > 4, and the corresponding RG 
equations also have an IR attractive fixed point in the 
physical range of the parameters. This establishes the 
scaling relations (5), (7), (8) with the new set of dimen- 
sions: 

A = 1 - 277/3 (exact), A e = A c = 2 - 277/3, 

A r = -2 + 2r?/3, A H = 2 + 277/3 (exact). (9) 

The dimensions A^/y are found exactly (there are no 
corrections of order rj 2 and higher) due to the Galilean 
invariance of the problem. In principle, the other dimen- 
sions are less universal than their analogues for the local 
case: besides the exponent 77, they can depend on d and 
a from (3). Our calculation has shown, however, that 
this dependence can occur only in the order 0(r/ 2 ). 

For d < 4, the model (4) with the nonlocal noise cor- 
relator ceases to be renormalizable: a new counterterm 
(</>') 2 is generated. A similar problem is well known in 
the RG approach to the stochastic Navier-Stokes equa- 
tion for purely incompressible fluid, where it occurs at 



d = 2; see e.g. Sec. 3.10 of Ref. [26]. In order to ap- 
ply the RG to this case, one has to extend the model 
by adding such term to the action from the very begin- 
ning, that is, one has to study the model with the mixed 
correlator 

D{k) =D Q + D' fc 4 - d - 2 " (10) 

(similar to that discussed in [13,14,17]). The extended 
model appears renormalizable, and its fixed points can 
be studied within the double expansion in two parame- 
ters, 77 and £ = 4 — d. 

The calculation in the first-order of such expansion 
shows that the extended model has two nontrivial fixed 
points. The first of them is IR attractive for e > 0, 
77 < 3e/4 and corresponds to the "local" regime with 
the dimensions (6), while the second is IR attractive for 
77 > 0, 77 > 3e/4 and corresponds to the dimensions (9). 
For 77 < 0, e < the only IR attractive point is trivial; it 
corresponds to a free (non-interacting) field theory. The 
regions of stability of the fixed points of the extended 
model in the e-77 plane are shown in Fig. 1. 

The main conclusions of our analysis are as follows. 
We have investigated a system of stochastic differential 
equations for the velocity and density of a self-gravitating 
matter, established two types of large-scale scaling be- 
havior (local and nonlocal ones), identified their regions 
of stability and calculated the scaling dimensions in the 
one-loop approximation (i.e., to first order of the corre- 
sponding e expansions). 

From the qualitative point of view, our analysis shows 
that nonequilibrium stochastic systems of the type (2) 
can have IR attractive fixed points in the physical range 
of parameters, and the corresponding scaling regimes can 
be treated systematically, within appropriate e expan- 
sions. What is more, such models can have several fixed 
points with different sets of dimensions, and the system 
undergoes the crossover in its large-scale behavior when 
its parameters (exponents in the forcing) change. 

It is worth noting that in model (2), the density and 
velocity fields have independent scaling dimensions, a fea- 
ture which is lost if the full set of equations is reduced 
to a single equation for only one independent field. Our 
results also suggest that rotational (non-potential) com- 
ponents of the velocity field do not decouple in those 
regimes and should be taken into account in the analysis 
of the large-scale behavior. 

Admittedly, the one-loop answers for the exponents are 
markedly larger than the latest experimental estimates 
for the exponent 7 (identified with 2Ag). In particular, 
for the local regime and d = 3 one obtains 7 « 3. For the 
nonlocal regime and arbitrary spatial dimension, 7 varies 
from 4 to 2 when the exponent 77 varies within its natu- 
ral range < 77 < 3/2 (for 77 > 3/2, the dimensions A^ >T 
become negative). This can be a hint that the simpli- 
fied model (2) does not include all physical interactions 
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relevant for the large-scale behavior (it is worth mention- 
ing here that the model (1) concerns directly dark mat- 
ter, while the galaxy-galaxy correlation functions concern 
visible matter). 

One can also expect that the simplest one-loop approx- 
imations (6), (9) overestimate the value of the scaling 
dimensions, and the second-order and higher corrections 
will improve the agreement, as indeed happens in the RG 
theory of fully developed turbulence; see [27]. Finally, it 
is possible that the scaling functions F in representations 
(5) are very singular in their arguments, which can lead 
to imaginary shift of the genuine exponent or to devia- 
tion from a plain power-law behavior, in agreement with 
some recent data [28]. 

In order to investigate these issues, one should go be- 
yond the simplest one-loop approximations and augment 
the plain RG equations by more advanced tools (renor- 
malization of composite operators, operator-product ex- 
pansion and so on), in analogy with the RG theory of 
fully developed turbulence; see e.g. [26]. This work is left 
for the future. 
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rj = 3e/4 




FIG. 1. Phase diagram of the model (4), (10), in the e-rj plane: the local regime is realized for e > 0, rj < 3e/4, the nonlocal 
one for rj > 0, ij > 3e/4 and the trivial one for rj, e < 0. 
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